Densest columnar structures of hard spheres from sequential deposition 
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The rich variety of densest columnar structures of identical hard spheres inside a cylinder can 
surprisingly be constructed from a simple and computationally fast sequential deposition of cylinder- 
touching spheres, if the cylinder-to-sphere diameter ratio is D £ [1, 2.7013]. This provides a direction 
for theoretically deriving all these densest structures and for constructing such densest packings with 
nano-, micro-, colloidal or charged particles, which all self-assemble like hard spheres. 
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I. INTRODUCTION 

The problem of finding the densest configurations of 
identical hard spheres inside a cylinder is related to many 
different aspects in science and engineering: Similar dens- 
est structures (e.g. zigzags and helices) have in recent 
years been observed for the self-assembly of nano-spheres 
inside carbon nanotubes [l|, charged particles on cylin- 
drical surfaces micro-spheres inside micro-channels 
[3] and colloidal particles inside cylindrical channels 3, 
and also for the formation of colloidal crystal wires a|) 
the reason being that in all such cases the particles close- 
pack themselves like hard spheres. Relevant problems of 
interest include (i) finding a theoretical derivation of the 
corresponding densest packings of hard spheres and (ii) 
finding a general methodology for constructing them in 
practice. 

For this problem, a wide range of chiral and achiral 
densest structures for a diameter ratio oi D = cylin- 
der diameter / sphere diameter G [1,1 + l/sm(7r/5)] 
have previously been found by simulated annealing [6|- 
|1| and by other computational means [1, [l3|, where 
1 -t- l/s«n(7r/5) — 2.7013 is a critical value of D below 
which the densest configurations are entirely made up 
of cylinder-touching spheres. In particular, the densest 
structures for D G [2,2.7013], which have their entire 
structural information contained on the cylinder surface, 
can be described by a phyllotactic scheme Q that is rel- 
evant to the study of plant morphology. 

In this article, it is shown that all the densest struc- 
tures for D G [1,2.7013] can surprisingly be constructed 
from a simple and computationally fast sequential depo- 
sition of cylinder-touching spheres which maximizes the 
number of spheres per unit length everywhere along the 
cylinder axis (typically less than 15 minutes of compu- 
tational time for each value of D; much faster than sim- 
ulated annealing which sometimes can take up a whole 
week 23 )• demonstrated in Fig. [T]for D G (1,2], the 
first four phases constructed from such a deposition pro- 



FIG. 1: The densest structures and their phyllotactic dia- 
grams for various cases of D € (1,2]: (a) the achiral planar 
zigzag at D = 1.800, (b) the chiral single helix at D = 1.940, 
(c) the chiral double helix at D = 1.992, and (d) the achi- 
ral doublets, which have two spheres located on one z posi- 
tion, at _D = 2.000. For each case, the required periodicity 
|V| — {D ~ 1)-K along s is indicated by a solid arrow. 



cedure are in complete agreement with the planar zigzag, 
single helix, double helix and achiral doublets predicted 
from simulated annealing 0, HI . The findings show that 
all these densest packings can be derived theoretically by 
means of maximizing the abovementioned number den- 
sity of spheres, and establish sequential deposition as a 
possible methodology for constructing such densest pack- 
ings with the aforementioned physical systems. 

This 3-D problem of packing identical spheres onto the 
inner surface of a cylinder is equivalent to a problem of 
packing objects in 2-D space, and is related to the 2-D 
problem of packing circular discs: In cylindrical polar 
coordinates (p, (f>, z) with a radial position p = (D — l)/2 
and an arc length s = (D — l)(/i/2, where the length 
quantities p, z and s are all expressed in units of the 
sphere diameter, the equation for two such spheres in 
contact is given by 
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FIG. 2: Plot of the vertical separation Az against the arc 
length As for spheres in contact at various values of the di- 
ameter ratio D. The regimes D £ [1,2) and D > 2 corre- 
spond to the absence and presence of a solution for Az = 0, 
respectively. As 7? — >■ oo, the relation between Az and As 
approaches (A^)^ -|- (As)^ = 1 for the 2-D packing of circular 
discs. 



which involves only two position variables (i.e. the differ- 
ences in z and (f> between the two spheres in contact), cor- 
responding to a 2-D problem, and can be derived easily 
using trigonometric identities. This implies the following 
(Fig. [2]): (i) For D S [1,2), the cylinder's surface is so 
curved that no two spheres can be accommodated on the 
same z position, i.e no solution exists for Az = 0. (ii) 
As Z? — J> oo, the cylindrical surface becomes flat and the 
equation for spheres in contact approaches that for the 
2-D packing of circular discs; In the densest configura- 
tion, the sphere centres form a perfect triangular lattice 
on the phyllotactic plot ^8] of z against s such that every 
sphere is in contact with 6 other spheres, (iii) For any 
finite value of D, the problem of finding an analytic so- 
lution to the densest configuration is complicated by the 
requirement of a periodicity |V| = [D — 1)tt along the s 
axis Q and by the relatively complex mathematical form 
of Eq. (m . Finding an analytic solution for cases of finite 
D thus remains a challenging theoretical problem. 

To gain insights into such an analytic solution, a semi- 
analytic approach is adopted here, based on the follow- 
ing conjecture: For any given value of D, if the packing 
of spheres is locally densest everywhere in the colum- 
nar structure, the overall structure would be at its glob- 
ally densest configuration. For a sufficiently long column 
of spheres, let AN/Az be the average number density 
of spheres (number of spheres per unit length) along z, 
where N are the indices of the spheres in ascending or- 
der of their z positions. This number density is related 
to the volume fraction, Vp = volume of spheres / volume 
of cylinder, as follows: 
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densest packing of spheres along z. It follows, from the 
abovementioned conjecture, that the densest columnar 
structure can be constructed from a sequential deposition 
of individual spheres onto the lowest possible z positions 
(determined across the entire 2tt range of (jj) at which 
they would come into contact with one or more other 
spheres. This is because such an iterative process would 
ensure the packing of spheres along z to be densest possi- 
ble at every z location. A template comprised of spheres, 
which covers the entire 27r range of (f>, is needed as a base 
for the initial deposition of spheres. For D G [1,2], Eq. 
(HI) for spheres in contact is valid for any value of A(j) 
so that, among other possible choices of a multi-sphere 
template, a single sphere is already enough to make up 
such a template. For D > 2, however, this equation is 
valid only for a limited range of Acj), which satisfies the 
condition cos Acj) G [1 — 2/(1? — 1)^, 1], so that a template 
comprised of more than one sphere is necessary. In gen- 
eral, for any given value of D, the columnar structures 
of spheres generated from the abovementioned deposition 
procedure depend on their underlying templates, and the 
densest structure can be obtained by varying across a di- 
versity of templates. 
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Thus the densest columnar structure, which by definition 
has the greatest volume fraction, also corresponds to the 



FIG. 3: Plot of the maximum volume fraction VF,max against 
the diameter ratio D, showing a complete overlap of data 
points between this deposition al gori thm and simulated an- 
nealing (data supplied by Mughal [ij]). The inset shows that 
the differences in VF,max between the two algorithms are of 
the order of 10"'^, i.e. less than 1% of percentage differences, 
and thus comparable to the numerical uncertainties in the 
simulations. The occasional peaks in |A1/F,maa:| are due to 
the relatively larger numerical uncertainties in the proximity 
of a structural transition. 



The different templates are generated through a depo- 
sition scheme in which the template structures depend 
only on the self-assigned relative positions (A(/))2,i = 
02 — '/'I € [0, tt] and (Az)2,i = Z2 — zi G [0, 1] between 
the two initially placed spheres: After the first sphere 
is placed at an arbitrary position on the plane of (s, z) 
coordinates, a particular value of (A0)2,i is chosen for 
the second sphere, where the value of (Az)2,i is given by 
the corresponding solution to Eq. ([Ij. If such a solu- 
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tion does not exist, the value of (Az)2,i is assigned sep- 
arately. Along the same deposition direction of increas- 
ing or decreasing (j>, which determines the chirality of a 
chiral structure, additional spheres are deposited if the 
existing template as made up of these two initially filled 
spheres does not cover the 27r range of (p. Inspired by the 
relevant triangular-lattice- like phyllotactic patterns Q, 
the positions of such additional spheres are chosen in the 
spirit to maintain a continuity in the relative positions 
A</)Ar,Ar_i = ipN - 't'N-i and Azn^n~i = zn - Zn-1 of 
successively filled spheres: For N > 2, the angular posi- 
tion of the A^*'' sphere is set to be {Acj))^,! = 4>n — 4>i — 
{N — l)(A0)2,i, and the sphere is deposited in a z po- 
sition at which it comes into contact with one or more 
other spheres. If such a z position does not exist, the 
sphere is placed at {Az)n.i = zn — zi = {N — l)(A2;)2,i. 
Deposition of such additional spheres stops whenever the 
template created has covered the 27r range of cj). The it- 
erative process of sequential deposition of spheres then 
starts from the angular position of the last deposited 
sphere in the process of template creation: In each iter- 
ative step, a sphere is deposited onto the lowest possible 
z position at which it would come into contact with one 
or more other spheres. Starting from the angular posi- 
tion of the previously filled sphere, this lowest z position 
is determined via a 27r scan along the same deposition 
direction of increasing or decreasing (jj as that for the 
process of template creation. Should there be degener- 
acy where this lowest z position occurs at more than 
one angular position, the first such angular position en- 
countered in the 27r scan is chosen. After a sufficiently 
long (e.g. 10 times of the sphere diameter) column of 
spheres has been generated, the average number density 
AN/ Az is evaluated via a linear fit of N against z, and 
the volume fraction is computed using Eq. The 
abovementioned procedure is repeated for a continuous 
range of (A0)2,i from to tt, and for a continuous range 
of (Az)2,i from to 1 whenever the value of (Az)2,i 
has to be assigned separately. As such, a quantitative 
agreement in the maximum volume fraction Vp^max be- 
tween simulated annealing and this deposition algorithm 
has been obtained for the whole range of D G [1, 2.7013], 
which is demonstrated for D e [1.750,2.7013] in Fig. 11 
(The numerical resolution of D is set to be 0.001, except 
for a better resolution near D = 2.7013). 

A complete agreement in the predicted densest struc- 
tures has also been obtained for D S [1,2.7013] between 
the two algorithms. For D e (1.000,2.000] the follow- 
ing sequence of densest structures has previously been 
predicted from simulated annealing 0, Q : (i) the achi- 
ral planar zigzag for D G (1.000,1.866], (ii) the chiral 
single helix for D e (1.866, 1.990], (iii) the chiral double 
helix for D e (1.990,2.000), and (iv) the achiral dou- 
blets, which have two spheres located on one z position, 
at D = 2.000. Here, the deposition algorithm yields the 
same sequence of densest structures for the same regimes 
of D. A gallery of the corresponding phyllotactic dia- 
grams is shown in Fig. [TJ In creating such diagrams, the 
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FIG. 4: (a) Sequence of structures as D increases from 2.040 
to 2.155: 321 symmetric at D = 2.040 Type 2 line slip, 
e.g. at D = 2.010 Type 3 line slip, e.g. at D = 2.142 -> 
330 symmetric at D = 2.155. (b) Two examples of periodicity 
breaking at _D = 2.350 and D — 2.620; in each case, the region 
that does not belong to the periodic structure is highlighted 
with a rectangle. 
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FIG. 5: Simultaneous deposition of m = 1,2,3 and 4 spheres 
at regular spacings of (j> in each deposition step. All the 
metastable curves for lower-density packings correspond to 
cases of u > 1. The symmetric structures are indicated along 
with their phyllotactic notations. 



boundaries of the equivalent close-packed objects in 2- 
D space are plotted from an equation which comes from 
the following transformation of Eq. ([2|): Az — )■ 2A2; 
and As 2As. This can be easily understood from 
the limiting case oi D — > 00: Applying such a trans- 
formation to (Az)^ -|- (As)^ = 1 for the 2-D packing of 
circular discs correctly yields the equation for a circle of 
radius = 1/2. For D e [2.000,2.7013], the phyllotactic 
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patterns from simulated annealing generally consist of 
regions with triads of interconnected spheres, which can 
be described in terms of a distorted triangular lattice [1] , 
as well as line defects that arise due to the required pe- 
riodicity |V| = (Z) — l)7r along s. These are referred to 
as line-slip structures At certain values of |V| or 
equivalently D, line defects are absent so that on a phyl- 
lotactic diagram every sphere appears to be in contact 
with 6 other spheres. Such defect-free structures, which 
generally correspond to the density maxima in Fig. ([3]), 
are described as symmetric Q. In such cases, the pe- 
riodicity vector is given by V = cibi -I- C2b2 where bi 
and b2 are any two of the three basis vectors. The coef- 
ficients ci > and C2 > are by definition chosen from 
the indices in the phyllotactic notation Imn where 
I = m + n. An example of such symmetric structures 
is the 220 structure in Fig. [Hd). For this regime of D, 
the densest structures predicted from the deposition al- 
gorithm are in complete agreement with the line-slip and 
symmetric structures predicted from simulated anneal- 
ing. For example, as shown in Fig. IUJa), the deposition 
algorithm yields the same sequence of structures as D 
increases from 2.040 to 2.155 0, @- 321 symmetric at 
D = 2.040 type 2 line slip for D e (2.040,2.142) 
-J> type 3 line slip for D e (2.142,2.155) 330 sym- 
metric at D — 2.155, where the different types of line- 
slip structures are defined in Ref. In some cases a 
periodic line-slip structure comes into existence only af- 
ter a few spheres have been deposited; Two examples of 
such periodicity breaking are demonstrated in Fig. Hl^b) . 
Lower-density packings, which might also form a subject 



of interest, can be generated from a simultaneous deposi- 
tion of multiple spheres at each deposition step (Fig. [5]) . 
For D > 2.7013, the structures with cylinder-touching 
spheres become hollow so that the aforementioned agree- 
ment with the 3-D densest structures from simulated an- 
nealing no longer exists. At D — !■ oo, the cylinder's sur- 
face becomes flat and essentially a 2-D triangular lattice 
of spheres is obtained. A sequential stacking of such 2-D 
layers along their perpendicular direction in AB or ABC 
sequence yields respectively the densest 3-D packings of 
HCP or FCC (Vf = 0.74048, the difference of which 
with the peak values in Figs. [3] and [S] clearly refiects 
the boundary effects imposed by the cylinder). For low 
values of D the experimental observations of the planar 
zigzag [ll, O-B, the helical phases P, and the achi- 
ral doublets [3, |1] could be related to the following: 

(i) The particles self-assemble like hard spheres, in the 
sense that their shapes remain approximately spherical. 

(ii) During the process of structural formation, the par- 
ticles are driven by a uni-directional force (e.g. gravity) 
which helps to maximize the number density AN/Az. 
In the future a more general deposition scheme that also 
applies to the regime of > 2.7013 and to the packing of 
binary hard spheres [H, [T!] should be developed. Possi- 
bilities of applying the current algorithm for the synthesis 
of chiral molecules should also be explored. 
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